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A family of dynamic models 
for large-eddy simulation 

By D. Carati 1 , K. Jansen, AND T. Lund 


1. Motivation and objectives 

Since its first application, the dynamic procedure has been recognized as an effec- 
tive means to compute rather than prescribe the unknown coefficients that appear 
in a subgrid-scale model for Large- Eddy Simulation (LES). The dynamic procedure 
(Germano et al. 1991; Ghosal et at. 1995) is usually used to determine the non- 
dimensional coefficient in the Smagorinsky (1963) model. In reality the procedure 
is quite general and it is not limited to the Smagorinsky model by any theoretical 
or practical constraints. The purpose of this note is to consider a generalized family 
of dynamic eddy viscosity models that do not necessarily rely on the local equilib- 
rium assumption built into the Smagorinsky model. By invoking an inertial range 
assumption, it will be shown that the coefficients in the new models need not be non- 
dimensional. This additional degree of freedom allows the use of models that are 
scaled on traditionally unknown quantities such as the dissipation rate. In certain 
cases, the dynamic models with dimensional coefficients are simpler to implement, 
and allow for a 30% reduction in the number of required filtering operations. 

2. Accomplishments 

2.1 A new family of dynamic eddy viscosity models 

The LES equations are obtained from the Navier-Stokes equations by applying a 
filter, denoted by an overline, which is assumed to damp scales smaller than A. In 
the context of eddy viscosity models, the unknown subgrid-scale stress generated 
by this_operation, r tJ = vJH J — u t Uj, is assumed to be proportional to the strain 
tensor Sij = ( d,Tij + d } Ui)/2: 

T ij ~ ~2,V t S ij. (1) 

The eddy viscosity, v t , has dimensions L 2 /T, where L is length and T is time. The 
characteristic length in the problem is obviously L c = A. Following the Kolmogorov 
(1941) dimensional analysis, the characteristic time may be expressed as a function 
of the rate of energy transfer within the inertial range £: T c = (A 2 /£‘) 1 / 3 . The 
“Kolmogorov expression” for the eddy viscosity is thus: 

v t = c k £ 1 / 2 A 4 / 3 , (2) 
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where c k is a non-dimensional constant. The rate of energy transfer is usually 
not directly accessible in LES, and thus Smagorinsky proposed to identify the rate 
energy transfer within the inertial range with the subgrid- scale dissipation: 

€ * -r, } S, } = * e |5| 2 , (3) 

where |5| 2 = 2S,jS lJ . When integrated over the volume, the above relation be- 
comes a good approximation since nearly all the dissipation will be carried by the 
subgrid-scale model when the cutoff is in the inertial range. In the Smagorinsky 
model, this equality is assumed to be valid at every point in space by invoking a 
local-equilibrium assumption between production and dissipation of energy. Insert- 
ing relation (3) into the Kolmogorov scaling for the eddy viscosity (2) gives the 
Smagorinsky model 


„ e = c,|5|A 2 , (4) 

where c 3 = c^ 2) is the non-dimensional Smagorinsky constant. In the Smagorinsky 
model, the time scale is seen to be |5| -1 . Thus, if local equilibrium is assumed, two 
expressions are available for the time scale in the eddy viscosity. By dimensional 
analysis, the eddy viscosity can depend on the ratio of these two time scales as well 
as on the fundamental scaling in Eq. (2). The most general model can therefore be 
written as 


v, = F (!^) £ l/5 A*/*, (5) 

where F is an arbitrary function. In particular, we may focus on a series represen- 
tation for F : 


|S| C ' £(l-0>/3 A (4+20)/ 3< (6) 

/=1 

Here are a sequence of numbers that define the exponents for the various terms 
in the series. They need not be integers. The parameters c/ are non-dimensional 
coefficients. As important special cases, note that n — l,Ci = 0 leads to the 
Kolmogorov scaling with c\ = c*, whereas n = l,Ci = 1 leads to the Smagorinsky 
model with C\ = c s . 

While Eq. (6) is rather general, it has the apparent drawback that the unknown 
dissipation rate, £, appears as a model parameter for ^ 1. Historically this 
defect has effectively excluded all models encompassed by Eq. (6) except for the 
Smagorinsky model. The situation has changed with the introduction of the dy- 
namic procedure, however, and it is possible to use Eq. (6) generally if it is recast 
in a slightly different form. If we assume that the test and grid filters are in the 
inertial range, then the dissipation rate as well as each of the model coefficients, 
c/, should be the same at two filtering levels. The product of the dissipation rate 
(raised to some power) and a model coefficient should also be invariant with filtering 
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scale, and thus the dynamic procedure may be used to determine the dimensional 
parameters c/ — Thus when Eq. (6) is recast in terms of c/, we can 

make use of Eq. (1) and write the subgrid-scale models at the grid and test level as 


n 


t„ = S,„ 

l-l 

{la) 

Tij = -2^Q|f| c 'A (4+2C,)/3 

(76) 


1=1 


where A is the test-filter width and Sij is the test-filtered strain rate. When 
Eqs. (7a) and (7b) are substituted into the Germano identity (Germano et ai 1991), 
a set of integral equations for the c/ are obtained. Following Ghosal et ai (1995) 
we can reduce the integral equations to algebraic relations if we constrain the co- 
efficients to have no spatial variation over the directions in which the test filter is 
applied. The end result is 


(M lk )c k = -(LamW), (8) 

where the Leonard tensor is given by Lij = UiU 3 — uiu 3 . The I th model tensor is 
defined as 

m<' ) = -2 (a< 4+2 <'>/ 3 - A (4+2<,) / 3 |f |G . (9) 

The left hand side of Eq. (8) is a matrix of products of these tensors: Af/* = 
Finally, () denotes a spatial average taken over the directions in which 
the test filter is applied*. Note that when n ^ 1, a linear system must be solved 
in order to determine the dynamic model coefficients. When the pure Kolmogorov 
scaling (n = 1, (i = 0) is used, the dynamic estimation for the eddy viscosity reduces 
to: 


u e ~ — 


1 (LijSij) 

2(a 4 / 3 — 1) (SijSij)' 


( 10 ) 


where a = A/A. This relation was derived earlier by Wong & Lilly, (1994). This 
model has the advantage that knowledge of the Smagorinsky time scale \S\ is not 
required, and thus the model is independent of the local equilibrium assumption. 
The Kolmogorov model also has the practical advantage that fewer filtering opera- 
tions are required as compared with the Smagorinsky model. This is true since 

the term |S|Sij does not appear in the Kolmogorov model. Finally, it should 


* In practice averaging is usually not performed in inhomogeneous directions even if these are 
included in the test filter. This inconsistency introduces an error that has been found to have a 
negligible impact on the simulation results (Ghosal et a/., 1995). 




Time, (Ui)/M 

FIGURE 1. Decay of resolved turbulent kinetic energy. : Dynamic Smagorin- 

sky model; : Dynamic Kolmogorov model; • : filtered experimental data of 

Comte-Bellot and Corrsin (1971). U is the mean advection speed in the wind tunnel 
experiments, M is the spacing between the bars in the turbulence-generating grid, 
and 0.5 q 2 is the total turbulent kinetic energy at the first measurement station. 



Wavenumber, kL/(2n) 

FIGURE 2. Velocity spectra. : Dynamic Smagorinsky model; : 

Dynamic Kolmogorov model; ■ , * : experimental data of Comte-Bellot and Corrsin 
(1971) for Ut/M = 98 and 171 respectively. L = 10.8 M is the length of a side of 
the computational box. The other scaling parameters are defined in Fig. 1. 
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be mentioned that models mixing the Kolmogorov and the Smagorinsky scalings 
(n > 2, Cl = 0, (2 = 1) could be investigated for situations with poorly developed 
inertial ranges. Indeed, in that case both Kolmogorov and Smagorinsky time scales 
might play independent roles and the dynamic procedure could determine the rel- 
ative weighting of these two scalings. 

2.2 Numerical tests 

As a first step in evaluating the new class of models, the Kolmogorov model 
(Eq. (10)) is tested in simulations of decaying isotropic turbulence. The simulations 
target the experimental measurements of Comte-Bellot and Corrsin (1971) and are 
performed with a pseudo-spectral code (Rogallo, 1981) using 32 3 mesh points. The 
equation for the model coefficient is averaged over the volume so that the coefficient 
is a function of time only. The simulations are initialized so that the 3-D energy 
spectrum agrees with the experimental data (up to the mesh wavenumber) at the 
first measuring station. The initial field is obtained by simulating the decay from an 
earlier time where the velocity phases are set at random. By iteratively adjusting 
the energy spectrum at the earlier time, it is possible to construct a field that has 
the desired energy spectrum as well as realistic phase information. The objective 
of the simulation is to predict the energy decay rate and the 3-D spectrum at the 
two subsequent experimental measurement stations. 

Figure 1 shows the kinetic energy decay history for the dynamic Kolmogorov and 
Smagorinsky models. There is little difference between the results of the two mod- 
els and both agree quite well with the experimental data. Near the starting point, 
the Kolmogorov model is seen to be slightly less dissipative than the Smagorinsky 
model. This could have to do with the fact that the initial field is generated with 
the Smagorinsky model and thus a transient is introduced when the model is sud- 
denly switched to the Kolmogorov scaling. Three-dimensional velocity spectra are 
shown in Fig. 2. Again there is very little difference between the two models. The 
spectra are seen to be slightly less damped at high wavenumbers in the case of the 
Kolmogorov model. This difference actually makes the Kolmogorov model agree 
slightly better with the experimental data at the final measurement station. 

The results of these tests suggest that the dynamic Kolmogorov model may work 
just as well as the Smagorinsky model. This is significant since comparable accuracy 
can be expected with 30% fewer filtering operations. The fact that the Kolmogorov 
scaling works also suggests that other terms in Eq. (6) may be useful as well. 

3. Future plans 

The Kolmogorov model will be tested next in turbulent channel flow. If is proves 
successful there it will be incorporated in the CTR complex geometry codes. Once 
these results are interpreted, we will study models that include more terms with 
the obvious first choice being a blend of Smagorinsky and Kolmogorov scaling (n = 

2, Ci = 0, C 2 = 1). 
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